0.1 Import Data

if (!require("pacman")) install.packages("pacman")
pacman::p_load(tidyverse, broom, skimr, devtools, ggpubr, mgcv, extrafont, mgcViz, here)

theme_set(theme_pubclean(base_size = 14))

# load the function to print GAM figures
source(here("R", "p_gam.R"))
# import
lesion <- read_csv(here("cache", "summarised_lesion_data.csv"))
weather <- read_csv(here("cache", "weather_summary.csv"))

dat <- left_join(lesion, weather, by = c("site" = "Location", "rep" = "Rep"))

0.2 Fit GAMs

For reproducibility purposes, use set.seed().

set.seed(27)

0.2.1 mod1 - s(Distance)

mod1 <-
   gam(
      mean_pot_count ~ s(distance, k = 5),
      data = dat
   )

summary(mod1)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08024    0.04751   22.74   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df    F p-value    
## s(distance) 3.926  3.996 78.4  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.482   Deviance explained = 48.8%
## GCV = 0.76522  Scale est. = 0.75394   n = 334
print(p_gam(x = getViz(mod1)) +
         ggtitle("s(Distance)"),
      pages = 1)

0.2.2 mod2 - s(Distance) + Precipitation

mod2 <-
   gam(
      mean_pot_count ~ sum_rain+ s(distance, k = 5),
      data = dat
   )

summary(mod2)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ sum_rain + s(distance, k = 5)
## 
## Parametric coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.124951   0.055685  20.202   <2e-16 ***
## sum_rain    -0.010853   0.007088  -1.531    0.127    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F p-value    
## s(distance) 3.926  3.996 78.71  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.484   Deviance explained = 49.1%
## GCV = 0.76443  Scale est. = 0.75086   n = 334
print(p_gam(x = getViz(mod2)) +
         ggtitle("s(Distance) + Precipitation"),
      pages = 1)

0.2.3 mod3 - s(Distance) + Windspeed

mod3 <-
   gam(mean_pot_count ~ mws + s(distance, k = 5),
       data = dat)

summary(mod3)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ mws + s(distance, k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.64401    0.11825   5.446 1.01e-07 ***
## mws          0.12273    0.03059   4.012 7.47e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F p-value    
## s(distance) 3.929  3.996 81.99  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.504   Deviance explained = 51.2%
## GCV = 0.73389  Scale est. = 0.72086   n = 334
print(p_gam(x = getViz(mod3)) +
         ggtitle("s(Distance) + Windspeed"),
      pages = 1)

0.2.4 mod4 - s(Distance) + Windspeed + Precipitation

mod4 <-
   gam(mean_pot_count ~ sum_rain + mws + s(distance, k = 5),
       data = dat)

summary(mod4)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ sum_rain + mws + s(distance, k = 5)
## 
## Parametric coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.674729   0.118523   5.693 2.78e-08 ***
## sum_rain    -0.014719   0.006968  -2.112   0.0354 *  
## mws          0.131154   0.030694   4.273 2.53e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##              edf Ref.df     F p-value    
## s(distance) 3.93  3.996 82.84  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.509   Deviance explained = 51.8%
## GCV = 0.72844  Scale est. = 0.71333   n = 334
print(p_gam(x = getViz(mod4)) +
         ggtitle("s(Distance) + Windspeed + Precipitation"),
      pages = 1)

0.2.5 mod5 - s(Distance + Windspeed) + Precipitation

mod5 <-
   gam(
      mean_pot_count ~ sum_rain + s(distance + mws, k = 5),
      data = dat
   )
## Warning in term[i] <- attr(terms(reformulate(term[i])), "term.labels"): number
## of items to replace is not a multiple of replacement length
summary(mod5)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ sum_rain + s(distance + mws, k = 5)
## 
## Parametric coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.124951   0.055685  20.202   <2e-16 ***
## sum_rain    -0.010853   0.007088  -1.531    0.127    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F p-value    
## s(distance) 3.926  3.996 78.71  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.484   Deviance explained = 49.1%
## GCV = 0.76443  Scale est. = 0.75086   n = 334
print(p_gam(x = getViz(mod5)) +
         ggtitle("s(Distance + Windspeed) + Precipitation"),
      pages = 1)

0.2.6 mod6 - s(Distance) + s(Windspeed) + Precipitation

mod6 <-
   gam(
      mean_pot_count ~ sum_rain + s(distance, k = 5) + s(mws, k = 5),
      data = dat
   )

summary(mod6)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ sum_rain + s(distance, k = 5) + s(mws, k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.18789    0.07152  16.608   <2e-16 ***
## sum_rain    -0.02613    0.01379  -1.895    0.059 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F  p-value    
## s(distance) 3.938  3.997 93.81  < 2e-16 ***
## s(mws)      3.932  3.997 15.97 3.05e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.566   Deviance explained = 57.8%
## GCV = 0.64971  Scale est. = 0.63051   n = 334
print(p_gam(x = getViz(mod6)) +
         ggtitle("s(Distance) + s(Windspeed) + Precipitation"),
      pages = 1)

0.2.7 mod7 - s(Distance) + s(Windspeed)

mod7 <-
   gam(
      mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5),
      data = dat
   )

summary(mod7)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08024    0.04366   24.74   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F  p-value    
## s(distance) 3.937  3.997 92.96  < 2e-16 ***
## s(mws)      3.917  3.995 16.04 6.03e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.562   Deviance explained = 57.3%
## GCV = 0.65392  Scale est. = 0.63659   n = 334
print(p_gam(x = getViz(mod7)) +
         ggtitle("s(Distance) + s(Windspeed)"),
      pages = 1)

0.2.8 mod8 - s(Distance) + s(Windspeed) + s(Precipitation)

mod8 <-
   gam(
      mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5) + s(sum_rain, k = 5),
      data = dat
   )

summary(mod8)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5) + s(sum_rain, 
##     k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08024    0.04345   24.86   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df      F p-value    
## s(distance) 3.938  3.997 93.802  <2e-16 ***
## s(mws)      2.977  3.051  6.713  0.0161 *  
## s(sum_rain) 1.931  1.942  1.979  0.1899    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.566   Deviance explained = 57.8%
## GCV = 0.64966  Scale est. = 0.63051   n = 334
print(p_gam(x = getViz(mod8)) +
         ggtitle("s(Distance) + s(Windspeed) + s(Precipitation)"),
      pages = 1)

0.2.9 mod9 - s(Distance) + s(Precipitation)

mod9 <-
   gam(
      mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5),
      data = dat
   )

summary(mod9)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08024    0.04594   23.51   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F  p-value    
## s(distance) 3.931  3.996 83.82  < 2e-16 ***
## s(sum_rain) 2.030  2.198 10.76 2.31e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.515   Deviance explained = 52.4%
## GCV = 0.71995  Scale est. = 0.70494   n = 334
print(p_gam(x = getViz(mod9)) +
         ggtitle("s(Distance) + s(Precipitation)"),
      pages = 1)

0.2.10 mod10 - s(Distance) +s(Precipitation) + Windspeed

mod10 <-
   gam(
      mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5) + mws,
      data = dat
   )

summary(mod10)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5) + mws
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  -2.8631     1.4637  -1.956   0.0513 . 
## mws           1.1095     0.4116   2.695   0.0074 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df     F  p-value    
## s(distance) 3.936  3.997 93.81  < 2e-16 ***
## s(sum_rain) 3.819  3.967 11.16 9.66e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.566   Deviance explained = 57.8%
## GCV = 0.64947  Scale est. = 0.6305    n = 334
print(p_gam(x = getViz(mod10)) +
         ggtitle("s(Distance) + s(Precipitation) + Windspeed"),
      pages = 1)

0.2.11 mod11.0 - s(Distance) + s(Windspeed) + s(Precipitation), family = tw()

This is the same as mod8 but using family = tw(), see ?family.mgcv for more on the families. The Tweedie distribution is used where the distribution has a positive mass at zero, but is continuous unlike the Poisson distribution that requires count data. The data visualisation shows clearly that the mean pot count data have this shape.

mod11.0 <-
   gam(
      mean_pot_count ~ s(distance, k = 5) + 
         s(mws, k = 5) + 
         s(sum_rain, k = 5),
      data = dat,
      family = tw()
   )

summary(mod11.0)
## 
## Family: Tweedie(p=1.044) 
## Link function: log 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5) + s(sum_rain, 
##     k = 5)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.22837    0.04099  -5.571 5.32e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df       F  p-value    
## s(distance) 3.497  3.855 123.761  < 2e-16 ***
## s(mws)      2.938  3.054  13.497 1.39e-05 ***
## s(sum_rain) 1.901  1.928   5.573   0.0132 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.673   Deviance explained = 61.2%
## -REML = 309.77  Scale est. = 0.36396   n = 334
print(p_gam(x = getViz(mod11.0)) +
   ggtitle("s(Distance) + s(Windspeed) + s(Precipitation), family = tw()"),
   pages = 1)

0.2.12 mod11.1 - s(Distance, bs = “ts”) + s(Windspeed, bs = “ts”) + s(Precipitation, bs = “ts”), family = tw()

mod11.1 <-
   gam(
      mean_pot_count ~ s(distance, k = 5, bs = "ts") + 
         s(mws, k = 5, bs = "ts") + 
         s(sum_rain, k = 5, bs = "ts"),
      data = dat,
      family = tw()
   )

summary(mod11.1)
## 
## Family: Tweedie(p=1.044) 
## Link function: log 
## 
## Formula:
## mean_pot_count ~ s(distance, k = 5, bs = "ts") + s(mws, k = 5, 
##     bs = "ts") + s(sum_rain, k = 5, bs = "ts")
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.22588    0.04094  -5.518 7.02e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##               edf Ref.df       F  p-value    
## s(distance) 3.250      4 117.921  < 2e-16 ***
## s(mws)      2.875      4  13.403  1.9e-12 ***
## s(sum_rain) 1.839      4   2.995 0.000496 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.671   Deviance explained =   61%
## -REML = 322.33  Scale est. = 0.36412   n = 334
print(
   p_gam(x = getViz(mod11.1)) +
      ggtitle(
         "s(Distance, bs = 'ts') + s(Windspeed, bs = 'ts')\n+ s(Precipitation, bs = 'ts'), family = tw()"
      ),
   pages = 1
)

This model, same structure as mod11.0, uses thin-plate splines to shrink the coefficients of the smooth to zero when possible.

0.3 Compare the Models

0.3.1 AIC, BIC

models <- list(mod1 = mod1,
               mod2 = mod2,
               mod3 = mod3,
               mod4 = mod4,
               mod5 = mod5,
               mod6 = mod6,
               mod7 = mod7,
               mod8 = mod8,
               mod9 = mod9,
               mod10 = mod10,
               mod11.0 = mod11.0,
               mod11.1 = mod11.1
               )
map_df(models, glance, .id = "model") %>%
   arrange(AIC)
## # A tibble: 12 x 7
##    model      df logLik   AIC   BIC deviance df.residual
##    <chr>   <dbl>  <dbl> <dbl> <dbl>    <dbl>       <dbl>
##  1 mod11.0  9.34  -320.  663.  708.     141.        325.
##  2 mod11.1  8.96  -322.  667.  712.     141.        325.
##  3 mod10    9.76  -392.  805.  846.     204.        324.
##  4 mod8     9.85  -392.  805.  847.     204.        324.
##  5 mod6     9.87  -392.  806.  847.     204.        324.
##  6 mod7     8.85  -394.  808.  845.     207.        325.
##  7 mod9     6.96  -412.  840.  870.     231.        327.
##  8 mod4     6.93  -414.  844.  874.     233.        327.
##  9 mod3     5.93  -416.  846.  873.     236.        328.
## 10 mod2     5.93  -423.  860.  886.     246.        328.
## 11 mod5     5.93  -423.  860.  886.     246.        328.
## 12 mod1     4.93  -424.  860.  883.     248.        329.

0.3.2 R2

enframe(c(
   mod1 = summary(mod1)$r.sq,
   mod2 = summary(mod2)$r.sq,
   mod3 = summary(mod3)$r.sq,
   mod4 = summary(mod4)$r.sq,
   mod5 = summary(mod5)$r.sq,
   mod6 = summary(mod6)$r.sq,
   mod7 = summary(mod7)$r.sq,
   mod8 = summary(mod8)$r.sq,
   mod9 = summary(mod9)$r.sq,
   mod10 = summary(mod10)$r.sq,
   mod11.0 = summary(mod11.0)$r.sq,
   mod11.1 = summary(mod11.1)$r.sq
)) %>%
   arrange(desc(value))
## # A tibble: 12 x 2
##    name    value
##    <chr>   <dbl>
##  1 mod11.0 0.673
##  2 mod11.1 0.671
##  3 mod10   0.566
##  4 mod6    0.566
##  5 mod8    0.566
##  6 mod7    0.562
##  7 mod9    0.515
##  8 mod4    0.509
##  9 mod3    0.504
## 10 mod2    0.484
## 11 mod5    0.484
## 12 mod1    0.482

0.3.3 ANOVA

anova(mod1,
      mod2,
      mod3,
      mod4,
      mod5,
      mod6,
      mod7,
      mod8,
      mod9,
      mod10,
      mod11.0,
      mod11.1,
      test = "F")
## Analysis of Deviance Table
## 
## Model  1: mean_pot_count ~ s(distance, k = 5)
## Model  2: mean_pot_count ~ sum_rain + s(distance, k = 5)
## Model  3: mean_pot_count ~ mws + s(distance, k = 5)
## Model  4: mean_pot_count ~ sum_rain + mws + s(distance, k = 5)
## Model  5: mean_pot_count ~ sum_rain + s(distance + mws, k = 5)
## Model  6: mean_pot_count ~ sum_rain + s(distance, k = 5) + s(mws, k = 5)
## Model  7: mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5)
## Model  8: mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5) + s(sum_rain, 
##     k = 5)
## Model  9: mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5)
## Model 10: mean_pot_count ~ s(distance, k = 5) + s(sum_rain, k = 5) + mws
## Model 11: mean_pot_count ~ s(distance, k = 5) + s(mws, k = 5) + s(sum_rain, 
##     k = 5)
## Model 12: mean_pot_count ~ s(distance, k = 5, bs = "ts") + s(mws, k = 5, 
##     bs = "ts") + s(sum_rain, k = 5, bs = "ts")
##    Resid. Df Resid. Dev          Df Deviance          F    Pr(>F)    
## 1     329.00     248.10                                              
## 2     328.00     246.34  1.00004467     1.76     4.8398  0.028515 *  
## 3     328.00     236.49  0.00032438     9.84 83384.0223 2.780e-11 ***
## 4     327.00     233.31  1.00009343     3.19     8.7512  0.003321 ** 
## 5     328.00     246.34 -1.00041782   -13.03    35.7856 5.797e-09 ***
## 6     324.01     204.37  3.99780494    41.97    28.8451 < 2.2e-16 ***
## 7     325.01     206.98 -1.00194276    -2.62     7.1733  0.007744 ** 
## 8     324.01     204.38  0.99830140     2.60     7.1576  0.007873 ** 
## 9     326.81     230.54 -2.79549912   -26.16    25.7123 3.516e-14 ***
## 10    324.04     204.44  2.76976222    26.11    25.8972 3.572e-14 ***
## 11    323.66     639.52  0.37229260  -435.08                         
## 12    323.68     644.19 -0.01578706    -4.66   811.7846 2.141e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.3.4 Check Best Model Fit

0.3.4.1 Mod11.0 - s(Distance) + s(Windspeed) + s(Precipitation), family = tw()

mod11.0_vis <- getViz(mod11.0)
check(mod11.0_vis,
      a.qq = list(method = "tnorm", 
                  a.cipoly = list(fill = "light blue")), 
      a.respoi = list(size = 0.5), 
      a.hist = list(bins = 10))
## 
## Method: REML   Optimizer: outer newton
## full convergence after 8 iterations.
## Gradient range [-0.0006677089,-7.453262e-09]
## (score 309.7676 & scale 0.3639647).
## Hessian positive definite, eigenvalue range [0.3939723,2978.454].
## Model rank =  13 / 13 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##               k'  edf k-index p-value   
## s(distance) 4.00 3.50    0.86   0.010 **
## s(mws)      4.00 2.94    0.89   0.045 * 
## s(sum_rain) 4.00 1.90    0.90   0.045 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.3.4.2 Mod11.0 - s(Distance, family = “ts”) + s(Windspeed, family = “ts”) + s(Precipitation, family = “ts”), family = tw()

mod11.1_vis <- getViz(mod11.1)
check(mod11.1_vis,
      a.qq = list(method = "tnorm", 
                  a.cipoly = list(fill = "light blue")), 
      a.respoi = list(size = 0.5), 
      a.hist = list(bins = 10))
## 
## Method: REML   Optimizer: outer newton
## full convergence after 8 iterations.
## Gradient range [-6.815834e-09,5.749888e-08]
## (score 322.3273 & scale 0.3641249).
## Hessian positive definite, eigenvalue range [0.6649335,2979.663].
## Model rank =  13 / 13 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##               k'  edf k-index p-value   
## s(distance) 4.00 3.25    0.86   0.010 **
## s(mws)      4.00 2.88    0.88   0.020 * 
## s(sum_rain) 4.00 1.84    0.89   0.025 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.4 Thoughts

This model, mod11.0, mean_pot_count ~ s(Distance) + s(Windspeed) + s(Precipitation) - family = tw(), is the best performing model. It cannot be used for predictions, but suitably describes the dispersal data we have on hand with the parameters used. More data would be desireable to increase the value of k as evidenced in the GAM checks.